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Abstract 

We review the problem of state reconstruction in classical and in quantum physics, which is 
rarely considered at the textbook level. We review a method for retrieving a classical state in 
phase space, similar to that used in medical imaging known as Computer Aided Tomography. We 
explain how this method can be taken over to quantum mechanics, where it leads to a description 
of the quantum state in terms of the Wigner function which, although may take on negative values, 
plays the role of the probability density in phase space in classical physics. We explain another 
approach to quantum state reconstruction based on the notion of Mutually Unbiased Bases, and 
indicate the relation between these two approaches. Both are for a continuous, infinite-dimensional 
Hilbert space. We then study the finite-dimensional case and show how the second method, based 
on Mutually Unbiased Bases, can be used for state reconstruction. 

PACS numbers: 03.65.Wj,03.67.Ac 
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I. INTRODUCTION 



The retrieval of the state of a physical system is an important problem in classical as well 
as in quantum physics, and yet it is a subject which is seldom discussed in textbooks. 

The state of a system in classical physics is described by a density in phase space, p(q,p), 
which could be determined (i.e., reconstructed) by the directly measurable conditional prob- 
ability density of the position q for a given momentum p, P(q\p), and the measurable prob- 
ability density of p, P(p), through the relation 

p(q,p)=P(q\p)P(p). (1.1) 

Thus we may envision the state as being specified by the set of measurable quantities P(q\p) 
and P(p). An alternative approach for determining a classical state involves measuring a 
linear combination of position and momentum (the constants a and b are introduced for the 
purpose of fixing dimensions) 

X e = aCq + bSp, C = cos6, S = sin6, (1.2) 

sometimes termed, for electromagnetic-field state measurements, "rotated quadratures" 
(Ref. p. 136). The probability for the new variable Xq for all values of 6 can then 
be used to reconstruct the phase-space density p(q,p) [2j. This procedure is similar to the 
familiar one employed in medical imaging for the reconstruction of a two-dimensional (2D) 
configuration density p(x,y), known as the Computer Aided Tomography (CAT) scan jj- 
4|: one simply replaces the two-dimensional configuration-space variables (x,y) of the CAT 
method by the two-dimensional phase-space variables (q,p). 

In quantum physics, a system may be prepared in a pure state described by a vector in 
Tilbert space, or, more generally, in a mixed state described by a density operator p (Ref. 
il, p. 204 and Ref. 6], p. 72). The problem of state retrieval involves the inverse inquiry, 
i.e., what are the measurable quantities whose values will suffice to determine the quantum 
state. Historically, this question may be traced back to the Pauli query j?J whether one 
can reconstruct the wave function, amplitude and phase, for a one-particle system, from the 
probability of its position, i.e., ^(rc)! 2 , and that for its momentum, |^(p)| 2 ; here ip(p) is the 
wave function in the momentum representation, the tilde indicating the Fourier transform. 
We now know that, in general, this is not possible: we need more information than these 
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two distributions. The literature on this subject, which is still of current interest, has grown 
enormously ever since. Here we have made a selection out of these approaches, with the 
idea of providing a link with the classical reconstruction scheme. 

The classical approach based on P(q\p), is, of course, untenable in quantum physics, where 
a fixed momentum precludes a well-defined position probability. A similar observation is 
applicable to the direct approach of measuring the joint probability of q and p. However, 
it is remarkable that the alternative method based on measuring Xg defined in phase space 
can be taken over to quantum mechanics (Ref. p. 143, Ref. js], p. 101). But then 
the question arises: how can that be, if there is no such thing as a joint probability density 
p(q,p) in quantum mechanics? It turns out that the answer one obtains by following this 
procedure is a function defined in phase space which, although is not a bona-fide probability 
density (it is real, but not-necessarily non-negative, and has sometimes been named a "quasi- 
probability" ) , contains all the information needed to compute any quantum mechanical 
expectation value we please, just as if we were given the complex wave function, or the density 
operator. This concept of quasi-probability was invented by Wigner 9] in the early days 
of Quantum Mechanics, with the purpose of finding the quantum-mechanical corrections 
to thermodynamic functions, and is known as the Wigner function. Thus retrieving the 
Wigner function using this tomographic method is a true quantum-state reconstruction, and 
to explain how this is achieved, and its relation with the classical tomographic approach, 
constitutes the main goal of the present paper. The main results for this approach are to be 
found in Eq. ( \3A\\ below for the classical case, and in Eq. ( 13. 16f) for the quantum-mechanical 
one. 

There is another concept which has been very useful in the task of reconstructing a quan- 
tum state. To give a trivial example, consider the eigenvectors of position and momentum: 
if the state vector of a system is an eigenstate of momentum, the system is equally likely 
to be found in any of the eigenstates of position. Pairs of bases with a similar property 
have been extensively studied 1(], 11] and are known as Mutually Unbiased Bases (MUB). 



It turns out that MUB constitute a powerful tool for state reconstruction, since it is possible 
to express the density operator that defines the state of the system in terms of a complete 



orthonormal set of operators 12Hl4j. We will explain the MUB approach to the problem 



of state reconstruction and show that the result [see Eq. (14.231) below] is consistent with 
that found with the method explained above, based on tomography in phase space and the 
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Wigner function. Even more important, we shall find that the two approaches correspond, 
essentially, to employing two ways of handling the same complete set of operators, thus 
providing a unified description of both methods. 

The paper is organized as follows. In the next section we review the CAT scan method, 
as employed for the reconstruction of a classical 2D density. In Section IIII Al we present a 
scheme for classical state reconstruction in phase space similar to the CAT method employed 
in configuration space. In Section IIII Bl we explain how the classical scheme can be taken 
over to quantum mechanics, and explain the role played by Wigner function. We then 
present in Sec. IIVI the alternative method for quantum state reconstruction based on the 
notion of MUB. So far, the discussion has been restricted to quantum systems described in 
a continuous, infinite-dimensional Hilbert space, because of our desire to make an analogy 
with classical physics. However, there have been many contributions to the problem of 
state reconstruction for quantum systems described in a finite-dimensional Hilbert space 
employing the notion of MUB. Although these systems do not have a classical counterpart, 
still they allow us to draw an illuminating parallel with the various concepts that have been 
introduced for a continuous Hilbert space. This fact motivates the brief discussion on the 
role of MUB for a finite-dimensional Hilbert space presented in Sec. |V] . Finally, we give our 
conclusions in Sec. IVII To avoid cluttering of the main text, we include some details of the 
mathematical derivations in a number of appendices. 

We wish to emphasize that the main goal for writing this paper is to give a pedagogical 
presentation of a subject which has been studied for many years and is still of current interest. 
With this motivation, we use a language that can be followed by a physics graduate student. 
We do hope that the analysis is in a form that allows its incorporation in a graduate Quantum 
Mechanics course. 



II. THE CLASSICAL RECONSTRUCTION SCHEME 

First we review briefly the method mentioned in the Introduction, the CAT scan, that 
is used for the reconstruction of a two-dimensional (2D) configuration density p(x,y). The 
mathematical procedure can be translated directly to retrieve a classical 2D phase-space 
density p(q, p) and, even more interesting for us, it can be taken over to quantum mechanics. 

In a 2D CAT-scan [2 - 4 ] a fine pencil beam of X-rays passes through a sample, shown as 
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the shaded area in Fig. [TJ along the "line of sight" defined by r • n = x' ; r is the position 
vector of a point on the "line of sight" and n a unit vector perpendicular to the line of 
sight, forming an angle 9 with the x-axis, so that r and n can be written asr = xi + ?/j, 
n = cos# i + sin 9 j, i and j being unit vectors along the x and y axes, respectively. Then 
the equation for the line of sight becomes 



The line of sight is offset by the amount x' Q from the rotated y' axis. The beam is attenuated 

y 



FIG. 1: In 2D Computer-Aided- Tomography (CAT), a beam of X-rays is passed through a sample, 
indicated by the shaded area, along the "line of sight" which is offset by the amount x' from the rotated 
y' axis. The unit vector n, perpendicular to the line of sight, forms an angle 9 with the x-axis. Knowing 
the response of the sample for all offsets x' and directions 9 we can reconstruct the original sample density 

p(x,y)- 

by scattering and absorption produced by the various parts of the sample encountered along 
the path. Assuming that the attenuation at (x, y) is proportional to the sample density 
p(x,y), the total attenuation will be proportional to 



where we have used an arbitrary offset value designated by x', and C and S have been defined 
in Eq. ( II. 2p . Now it is important to remark that knowing the response of the sample given 
by po(x') for all x' and directions 9, we can reconstruct the density p(x,y) of the sample. 



x' = x cos 9 + y sin 9. 



(2.1) 




line of sight 




(2.2) 
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The mathematics of this problem was actually developed by J. Radon at the beginning of 
the twentieth century [ij]] for the study of astronomical data. In fact, the function pe(x) of 
Eq. (12. 2p is known in the literature as the Radon transform of the density p(x, y). Thus the 
task is to invert the Radon transform to find the sample density. 

It is shown in Appendix |A] that the sample density p(x, y) can be expressed in terms of 
the response of the sample, po(x'), for all x' and directions 6 defined above (see Ref. [lj], pp. 
144), as 

pf Xty) = _j_ r dd v r dx> dMx ' )/dx ' , (2.3) 

2ir 2 Jq J_ 00 x' — (x cos 9 + y sin 6) ' 

where V stands for the Cauchy principal value of the integral. Indeed, Eq. ( 12. 3ft is the 
inverse Radon transform of po{x'). 

To gain some insight into the structure of the sample response pe(x'), it is illustrative to 
consider the particular case in which the sample density p(x,y) is isotropic, i.e., dependent 
only on the distance r = x 2 + y 2 from the origin and independent of the angle. If we write 
x and y in polar coordinates as x = rcos0, y = rsin0, Eq. (12. 2p for the response pe(x') 
takes the form 

Pe(x') = dcj) I dr r 5(x' — r cos(0 — 6))p(r), 
Jo Jo 

2tt poo 

dr r 5(x' — r cos 0)p(r), (2.4) 

showing that pe(x') is independent of 9 for the isotropic case. By direct substitution, one 
may also observe that in this case pe{%') is symmetric, i.e., pg(-x') = pe(x'). 

It may also be pointed out that the sample density p(x, y) must be a non-negative quantity, 
although this fact is not explicitly manifest in Eq. ( 12. 3ft . It is thus useful to verify this 
property in some particular example. For this purpose, we choose the isotropic case studied 
in the last paragraph. For example, at the origin of coordinates, x = y = 0, Eq. (12. 3 p gives 

pM= -^ v jy°j^, ( , 5) 

Since pe(x') is a symmetric function of x', dpg(x')/dx' is antisymmetric. Since the quantity 
—x' appearing in Eq. (12. 5p has precisely this same property, as illustrated in Fig. [21 the 
resulting density p(0, 0) at the origin of coordinates is positive. 
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FIG. 2: Schematic representation of the functions pg(x'), dpg{x') / dx' (solid lines) and —1/x' 
(dashed line) as functions of x' , needed in the text to show that the sample density /)(0,0) at 
the origin is positive. 



III. CLASSICAL-QUANTUM PHYSICS STATE-RECONSTRUCTION ANAL- 
OGY 

A. Classical state reconstruction 

A state in classical statistical physics is determined by a probability density in phase 
space. In this paper we shall always consider, for simplicity, one-particle systems with one 
degree of freedom. We write the probability density in phase space as p{q,p) which, for 
convenience in our comparison with Quantum Mechanics, will be normalized as 



After the discussion given in the previous section on CAT in 2D configuration space (x, y), 
it is clear that a similar method can be applied in 2D phase space (q,p): if we consider the 
linear combination of position and momentum given in Eq. (11.2)1 . the probabilities for the 
new variable Xg for all values of 6 can then be used to reconstruct p(q,p) 3|. 



Before proceeding, we indicate our choice for the constants a and b which were introduced 
to fix dimensions. We choose a = 1/qo, b = 1/po, where qo and po represent any convenient 
scales for position and momentum. Subsequently renaming the dimensionless quantities 




(3.1) 
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q/qo and p/po again as q and p, respectively, the transformation of Eq. (I1.2p reads 

X e = Cq + Sp. (3.2) 

We go back to the probability density of the variable X 9 . If we designate it as pe(x'), 
where x' represents an arbitrary value of X e , we have, just as in Eq. (12. 2p 

pe(x') = j j 5(x' -Cq- Sp) p(q,p) ^ . (3.3) 

The goal is to find p{q,p) in terms of pe{x') by inverting Eq. (13. 3p . Proceeding as in the 
previous section and Appendix |AJ we find the equivalent of Eq. (12. 3p as 

^ )= -ir^r^ , f^f . .. . 0.4) 

71 Jo J-oo x {Q cos v + p sin v) 

B. Quantum state reconstruction 

As mentioned in the Introduction, the above method based on measuring Xg defined in 



phase space can be taken over to quantum mechanics (Ref. p. 143). This leads to a 
quasi-probability density defined in phase space known as the Wigner function. 

In what follows we shall take units in which h = 1. Consider an arbitrary Hermitean 



operator A. We define its Wigner transform as 



y 



W A (q,p) = J e~™{q+»- A\q- V -)dy. (3.5) 

For the case where the operator A is the density operator p defining the state of the system, 
we speak of the Wigner function of the state, which has the normalization property 

/ / Wfap)^ = 1 • (3-6) 

similar to the normalization of Eq. (13. ip adopted for the classical distribution. 

It is well known [l| that Wigner function for a state may be negative in some parts of 
phase space. Thus it does not qualify as a true probability density and is referred to as a 
quasi-probability density. An illustration of the fact that it plays in quantum mechanics a 
role analogous to that played by the classical probability density p(q,p) is the similarity of 
Eqs. (I3.15P and (13 . 1 6[) given below with Eqs. ( 13 .3p and (13. 4p . respectively. 

An important property of Wigner function, obtained from the definition (13. 5p . is jl, 8| 

Tr{AB)= ( f W A (q,p)Wz(q,p)^ , (3.7) 



for any two operators A and B. This implies that the trace of the product of two operators 
in Hilbert space can be evaluated as an integral in phase space of the corresponding Wigner 
transforms. The normalization of Eq. ( 13. 6 j) is consistent with the property (13. 7p . taking 
A = p and B — 1. The statistical expectation value of an observable A, obtained by using 
Eq. (13. 7p . can be expressed as 

(A) = Tr(pA) = J J W fap)W x {q,p)^ , (3.8) 

i.e., as an integral in phase space of Wigner function for the state times Wigner transform of 
the observable. With these results, Wigner function of the state and the Wigner transform 
of observables can be employed to "do QM in phase space" . 

It is also a simple exercise to show that the above definition of the Wigner function of 
the state p is equivalent to the inverse Fourier transform of the characteristic function of the 
density operator (Ref. js|, Eqs. (3.12), (3.16)), i.e., 

Wp(q,p) = — [ I W(u,v)e i{uq+vp) dudv , (3.9a) 



py ~ lir ' 2ix j 

W(u,v) = Tr[p e ~ i{u{i+vi>) } . (3.9b) 



Now consider the observable 

X e = Cq + Sp, (3.10) 

which is the QM counterpart of the classical quantity of Eq. (13. 2h . This observable satisfies 
the eigenvalue equation 

X e \x';6) = x'\x';6} , (3.11) 

where x' denotes an eigenvalue and \x'; 9) the corresponding eigenvector. 

Our program is as follows. If the system is prepared in the state defined by the density 
operator p, we first consider the probability density p^ M (x') that a measurement of the 
observable Xg will give the value x'\ this probability density will be initially expressed in 
terms of p, Eq. (13.121) below, and then in terms of the Wigner function Wp(q,p) in phase 
space, Eq. (13.151) below. The final goal is to "invert" this relation and show that we can 
retrieve the Wigner function in terms of p® M \x'). 

The probability density pjj M {x') is given by the standard QM expression 

pT\x') = Tr(p P^,) , where F e x , = \x'; 6) (x'; 9\. (3.12) 
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Making use of Eq. (13. 8 j) , we write 

p ° M{x ' ] = J I w ^ w K^P) d ~^ ■ ( 3 - 13 ) 

In this expression, W ¥ e (q, p) is the Wigner transform of the projector P^,, which is calculated 
in Appendix IB1 with the result 



i (q,p) = 5(x'-(Cq + Sp)). (3.14) 

X' 

Then Eq. (13.131) takes the form 

pf \ x >) = J J Wp {q,p) 5{x' - (Cq + Sp)) ^ • (3.15) 

This last equation is the QM counterpart of Eq. (13 .3p for the classical probability density 
pe{x'). It shows explicitly that what plays the role of the classical probability density p(q,p) 
in phase space is now the quasi-probability density given by the Wigner function Wp(q,p). 
Thus in order to invert Eq. ( 13 . 1 5|) we just copy the result in Eq. (13. 4p and write Wp(q,p) in 
terms of p® M (x') as 

W p (q,p) = - 1 - fdBV r dx> , d f M ^l dx [ . (3.16) 
Jo J-oo x' - (qcosO +psm9) 

This equation allows reconstructing the QM state, in the sense that from the observable 
probability density pjj M (x') the Wigner function of the density operator can be extracted; 
its knowldedge, in turn, is equivalent to that of the state itself. 

This completes our analysis that shows a close analogy between the classical and quantum 
state reconstruction: both require the use of the inverse Radon transform. We now turn to 
an alternative quantum-state reconstruction scheme which does not require the use of the 
Radon transform. 



IV. MUTUALLY UNBIASED BASES AND STATE RECONSTRUCTION 

Mutually unbiased bases (MUB) in concept were introduced by Schwinger [l(| in his 
studies of vectorial bases for Hilbert spaces that exhibit "maximal degree of incompatibility" . 
The eigenvectors of x and p, \x) and \p), respectively, are example of such bases. The 
information-theoretical oriented appellation "mutual unbiased bases" was introduced by 



Wootters 



111 
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Consider two complete and orthonormal vectorial bases, Bi, B2, whose vectors will be 
designated by \u;Bi) and \v\B 2 ), respectively. The two bases are said to be MUB if and 
only if, for B\ 7^ B 2 , 

I (it; B 1 \v;B 2 )\ 2 = K , Vu,v, (4.1) 

where K is a constant independent of u and v (see Ref. [3]). This property means that 
the absolute value of the scalar product of vectors from different bases is independent of the 
vectorial label within either basis. This implies that if a system is measured to be in one 
of the states, say \u;Bi), of Bi, it is equally likely to be found in any of the states \v; B 2 ) 
of any other basis B 2 , when B\ and B 2 are MUB. The value of K may depend on the bases 
£>i, B 2 , which indeed is the case for a continuous Hilbert space. For a Hilbert space with a 
finite dimensionality d, K = 1/d. 

The concept of MUB is found to be of interest in several fields. For instance, the ideas 
are useful in a variety of cryptographic protocols {3] and signal analysis 

In what follows we outline a scheme for state reconstruction based on MUB jl8] which is 
an alternative to the one presented in the previous section. 



A. Some properties of the operator Xg and its eigenstates 

We begin with a review of the properties of the operator Xg and its eigenstates \x',9), 
Eqs. (13301) . (EOT]) , and show that the bases {\xi,9i)}, {\x 2 ;6 2 )} (6>i ^ 9 2 , fixed) are MUB. 
We repeat the definition (13.101) of the operator Xg and introduce the new operator Pg as 

Xg = Cx + Sp, Pg = -Sx + Cp ; (4.2) 

Xg and Pg are canonically conjugate, i.e., [Xq^Pq] = i, just as the original operators x, p. 

As a first step we solve the eigenvalue equation (13.111) in the coordinate representation. 
In this representation we define the wave function 

ij x ,,g(x) = (x\x';9), (4.3) 

which satisfies the equation 

^x cos 9 - i sin 9 -^j ^ x >,e{ x ) = x> ^x'fiix). (4.4) 
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The solution of this equation is 

= F(x',9) e -^ {x2cosd - 2xx '\ (4.5) 

where F(x', 9) is an arbitrary function of x' and 9. It is shown in Appendix O that F(x', 9) 
can be completely determined, up to an arbitrary overall phase, by imposing on the states 
\x';9) the requirements [2^] 

(xi, 9\x 2 , 9) = 5(xi - x 2 ), (4.6a) 
(xi, 9 \X g \x 2 , 9) = x 2 8(x 1 -x 2 ); (xi, 9\P e \x 2 , 9) = -i 8'(x 1 - x 2 ), (4.6b) 
^x',e( x ) = (x\x',9) ->• 8(x - x'), as 9^-0; ^ x ',e=7r/2{x) = —7= ■ (4.6c) 

V Z7T 



Here, Eq. ( I4.6ap expresses the ortho-normalization (in the sense of the Dirac delta func- 
tion) of the states \x';9). Equation ( I4.6b|) requires that the matrix elements of the new 



canonically conjugate operators Xg and Pg with respect to the new states \x'; 9) be equal 
to the matrix elements of the old canonically conjugate operators x and p with respect to 
the old states \x), as demanded by a canonical transformation. The first Eq. fl4.6cl) requires 
that in the limit 9 — > the overlap between the new state \x', 9) and the old one \x) be a 
Dirac delta function. The second Eq. f)4.6cp requires that for 9 = ir/2, i.e., when X e=n / 2 is 
the momentum p, the wave function ip x ',e=n/2{%) be a plane wave with no extra phases. 
The final result for the wave function ipx',e{x), up to an overall constant phase, is 

iff sgn(sin0)-|l 

\/2ir\ sin 6*| 

Notice the symmetry of this expression under the interchange x x'. Since | — x', 9 + 71) = 
\x', 9), it suffices to consider state vectors in the range —00 < x' < 00 and —ir/2 < 9 < ir/2, 
other values of 9 repeating the eigenvectors in this range (see Ref. p. 144). 

We relate the new state \x', 9) to the old one \x') through a unitary transformation as 

\x',9) = U\9)\x'). (4.8) 

For the reader's convenience, we mention that the operator U used in this article coincides 



with the one designated by V in Ref. 19( , and that called W in Ref. 20] . Using Eqs. (14. 8 p 
and (14.71) we find, for the matrix elements of the unitary operator U'(9) in the old basis, 

i[fsgn(sin0)-f] 

(x\U j (9)\x') = , e -^[(x'^)cose-2xx']_ ( 4 _g) 

a/27t| sin#| 
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Using the unitary operator U (9) we write the eigenvalue equation (13. lip as 
U (9)XqU\9)\x') = x'\x'), implying x = U{9)XeW (9). Thus the operator x and, similarly, 
its canonically conjugate p transform as 

X e = U\6) x U (6), and P e = U\9) p U (9) . (4.10) 

The above unitary transformation is given by the operator 

U{9) = e~ m , (4.11) 

where n = a^a is the number operator, and a, a) are the annihilation and creation operators, 
respectively, given by a = -^(x + ip) , at = -^(x — ip) . Indeed, using the operator identity 
(Ref. 3, P- 339 ) 

e A Be~ A = B+[A,B] + [A, B}} + ■ ■ ■ , (4.12) 
we readily find that the operator (14. lip gives the transformation properties of a and a) 

U\6)aU(0) = e- l8 a, U\9)tfU{9) = e ie a ] , (4.13) 

that lead to the transformation properties of x and p, Eq. ( I4.10p [with Eq. ( 14.21) ] . It is 
also shown in Appendix |D] that the matrix elements of the operator exp (i9n) (the adjoint 
of (14. lip ) are identical to those of Eq. ( 14. 9 p found earlier. 

Finally, we verify that the bases and {\x2~-, ^2)} (with fixed 9\ 7^ 62) that we have 

been studying above are MUB. Using Eq. gjj]) and the relation U(9 2 )U ] (9 l ) = U ] (9 X - 9 2 ), 
which follows from ( 14.111) . we find 

\{x 2] 9 2 \x 1 -9 l )\ 2 = \(x 2 \U\9 l - 9 2 )\ Xl )\ 2 = 27r|g( ^ A)| , (4.14) 

where S(9i,9 2 ) = sin(6 l 1 — 9 2 ). The number \(x 2 \9 2 \xi\9 1 )\ 2 is thus independent of X\ and 
x 2 , so that, according to definition ( 14 -If) , the two bases are MUB. As an example, for 
9 = 7r/2, \x';9 = ~) is an eigenfunction of p with eigenvalue x', whose projection in the 
x representation is e lx ' x /\/2n [see the second Eq. (I4.6cl) ]. its absolute value squared being 
consistent with Eq. ( 14. 14[) . 

In Appendix [E] we present a simple way to derive the result of Eq. (14.141) , which is an 
application of the idea of doing QM in phase space using Wigner transforms, mentioned 
right below Eq. ( 13Tg|) . 
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B. State reconstruction based on MUB 



Now we show that the MUB introduced above can be used to perform a quantum- 
mechanical state reconstruction. We first introduce the set of operators 



Z(a,6) 



e ia± e ibp 



kab i(ax+b 



e 2 e 



-oo < a, b < +00, 



(4.15) 



where we have used the BCH identity, Eq. ( 1B4I) . These operators form a complete and 
orthogonal operator basis [l, 8|. They satisfy the orthogonality property 

da db 



Tr 



Z\a',b')Z{a,b) 



2tt 



8{a' -a)5{b' -b). 



(4.16) 



Thus we express the density operator as a linear combination of the operators Z(a, b) as 

da db 



p — c(a, fe)Z(a, &)- 



2tt 



c(a,b) = Tr pl)[a,b) 



(4.17) 



OO 

dr r 



-ir(Ci+5p)] ir(Cz+Sp) 



(4.18) 



In the above equations, a and 6 play the role of Cartesian coordinates. We go over to polar 
coordinates, defining a = r cos 9, b = r sin 9, so that Eq. (I4.17P takes the form 
A _ 1 

P 27T j 

[We use the abbeviations C and S from Eq. (11. Using similar arguments to those that 
led from Eq. (IA5j) to flA6j) . we rewrite the above equation as 

1 



Tr [p 



P 



dt \t 

27T ./q 



Tr [p 



(4.19) 



Since in the exponent of this last equation we have the operator Xg [see Eq. f)4.2p ]. the 
exponential can be written in its spectral representation [see Eq. (13. lip 



as 



Jt(Cx+Sp) 



JtXa 



e ltx 'F x , dx' 



where the projection operator P x /,0 is defined in Eq. (13. 1 2D . Similarly, 



Tr [p 



-it(Cx+Sp) 



-itx 



Tr(pF x ',9)dx' 



- itx 'pQ M (x')dx', 



(4.20) 

(4.21a) 
(4.21b) 



where we have used the definition of the QM probability density p® M (x 1 ), Eq. (I3.12p . Then 
Eq. P~T9|) for p, using Eqs. KM and <14.21bj> . becomes 

j- j dt\t\J d9 j j dx'dx"e- it{ - x '- x "^p^ M {x'% 



P 



— ' lim / d9 

27T e^0+J 




dx'dx"f € {x' - x")p% M {x')F x „, 



(4.22a) 
(4.22b) 
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In the last line we have performed the radial integral and used the definition f ]A8j) that 
was introduced in our earlier analysis, in the course of inverting the Radon transform. It 
is important to note that in the present context we have been able to express the density 
operator p directly in terms of the probability pg(x'), thanks to the expansion of p, Eq. 
f l4.19p . in terms of MUB, together with Eq. ( 14.21b[) which relates the trace on its left-hand 
side with p® M \x'). 

The Wigner function for the state p of Eq. f!4.22bj) is identical to the result found above 
in Eq. ( 13 . 1 6j) . which we reproduce here for completeness 



W~ p {q,p) 



7T 



- d9V 



dx' 



dpj M [x')ldx' 
x' — (q cos 9 + p sin 6) 



(4.23) 



This result can be proved as follows. Application of Eq. ( 13. 5 p -defining the Wigner function- 
to the density operator p, Eq. ( I4.22b|) . gives 



/»7r /> r /"OO 

Wp(q,p) = — lim / d9 dydx'dx" f t {x' 

2tt e ->o+ Jo J J 



^~ QM (x') 



xe 



-ipy 



P, 



y 

q 

y 2 



(4.24) 



We evaluate the matrix element of the projector P x ",0 by using its definition in Eq. (13. 12ft . 
the unitary transformation, Eq. ( 14. 8p . and its explicit expression, Eq. ( 14. 9L to find 



!J 



Q 



y 



(4.25) 



(4.26) 



2/ 27rsin# 

Substituting this result in Eq. (I4.24p and performing the integration over y we have 
Wp(q,p) = lim — d6 dx'dx"f e (x' -x")p$ M (x') 

e^0+ 2li J J J^ 

xS(x" — (gcos^ + psin.9)) 

= lim — d9 dx' f e (x' — (gcos 6* + psin8))p9 M (x) . 
e^o+ 2tx Jo J-oo 

The last line is 2ir times the right-hand side of Eq. (IA9[) . with x replaced by q and y by p. 
We thus take over the result of Eq. ( 12.31) . making these replacements and multiplying by 
2tt, and find Eq. (T4T231 . 

Finally, we calculate the matrix elements of the density operator (14.22b[) in the coordinate 
representation, (xi\p\x2), which is the counterpart in Hilbert space of Eq. (I4.23p . For the 
matrix elements of the projector IV',0 we find, just as in Eq. (I4.25p . 



X-i 



+g(x 1 -x 2 ){x 



II x 1 +x 2 



cos 9) 



2ir sin 9 



(4.27) 
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so that 

1 r n r r°° p^?( xi ~ ^x^"— xi ^ x2 cos8) 

( Xl \p\x 2 ) = — \im d6 dx'dx"f € (x'-x")p% M (x') — — — . 

(4.28) 

We compare this last equation with Eq. (IA9[) and use the result of Eq. (I2.3P to obtain 

1 r r r°° dr'dr" dn QM (r') p^tei^-^)^"-^ 2 - cos6 ) 
<*i|/W = — d6V / , P * X) - , (4.29) 

which shows explicitly how p® M (x'), which is a probability density, and hence a measurable 
quantity, can be used to find the matrix elements of the density operator. 

This completes our demonstration of the consistency of the two approaches to the prob- 
lem of quantum-mechanical state reconstruction that we have considered in this paper, for 
systems described in a continuous Hilbert space. On the one hand, the approach presented 
in the previous section based on tomography in phase space and the Wigner function and, 
on the other, the one given in the present section based on the expansion of the density 
operator in terms of operators defined via MUB. The difference in the strategies of these 
two approaches involves, essentially, two ways of handling the complete orthonormal oper- 
ators e iax e lb P; the Wigner function approach that led us to the results in Sec. IIIIBI can be 
regarded as using these operators to construct the Fourier transform of the density operator, 
as in Eq. (I3.9p . If, on the other hand, we consider their spectral representation, Eq. (I4.20p . 
we are led to the MUB approach of the present section. 



V. MUTUALLY UNBIASED BASES AND STATE RECONSTRUCTION IN A 
FINITE-DIMENSIONAL HILBERT SPACE 



Considerab 
Hilbert space 



e work has been devoted to the study of MUB in a finite, <i-dimensional 
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LH24J. In this paper we restrict our study to the case in which the dimen- 
sionality d is a prime number, for this case the number of MUB is exactly d + 1 [ll|, 



14 



22]. 



The finite- dimensional theory is intriguingly connected with sophisticated mathematical no- 
tions 
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25] that we do not consider here. 



In the finite, d- dimensional Hilbert space problem, where a Radon-like transform is not 
available, we shall follow a procedure which is analogous to that presented in the last section 
for a continuous, infinite-dimensional Hilbert space. 
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We first consider the <i-dimensional Hilbert space to be spanned by d distinct states 
|n), with n — 0, 1, • • • , (d — 1), which are subject to the periodic condition \n + d) = \n). 
These states are designated as the "computational basis" of the space. We shall follow 



Schwinger 10( and introduce the unitary operators X and Z, which play a role analogous to 
that of the position operator x and the momentum operator p of the continuous case. The 
Schwinger operators are defined by their action on the states of the computational basis by 
the equations 

Z\n) = u n \n), to = e 2m/d , (5.1a) 
X\n) = \n + l). (5.1b) 

These definitions lead to the commutation relation 

ZX = uXZ. (5.2) 

The two operators Z and X form a complete algebraic set, in that only a multiple of 

n 

the identity commutes with both [10(. As a consequence, any operator defined in our d- 
dimensional Hilbert space can be written as a function of Z and X. 

The d 2 - dimensional matrix space is spanned by the complete orthonormal d 2 operators 
X m Z l , with m, I = 0,1,. .{d — 1), so that any d x d matrix can be written as a linear 
combination of these d 2 operators. A familiar example is a 2-dimensional Hilbert space, 
where any 2x2 matrix can be written as a linear combination of the three Pauli matrices 
plus the unit matrix, which can also be written as a x , a z , u x a z and I. 

The operators X m Z l are orthonormal under the trace operation, 



Tr 



x m z l (x m 'z l " 



d 5 mm i5n>, (5.3) 



a relation which can be proved directly using the defining Eqs. ( 15. ip . Completeness follows 
from the set consisting of d 2 linearly independent operators. 

We shall need to do arithmetic operations on the numbers n = 0, 1, • • • , {d — 1) that label 
our states, assuming the periodic condition d = 0[mod d}. When d is a prime number, the 
operations of multiplication and division, modulo d, can be defined consistently {26 1. As a 
simple example, we find, for d = 3, that 1/2 [mod 3] = 2, since 2-2 = 4 = l[mod 3]. 

As a result, we may replace uniquely, up to a power of u, all the operators of the form 
X m Z\ with m ^ 0, by (XZ b ) m , with b = lmr x = 0, 1, ...{d - 1). We can readily prove that 

X m Z l = uj~ v ^ }1 \XZ h ) m , (5.4) 
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with I = mb[mod d]. ^From Eq. (15. 3p we find, for the new quantities (XZ b ) m , the orthogo- 
nality relation 



Tr 



i xz'T ( (xz 



'rb'\m' 



At 



b,b' Vm,m' i 



m, m! ^ 0. 



(5.5) 



The operators (XZ b ) m are d(d — 1) in number; these, together with the d operators Z l , 
with I = 0, 1, ...(d — 1) (this last set evidently contains the identity: Z° = Z d = I), form 
a complete orthonormal set of d 2 operators which is equivalent to the set X m Z l considered 
above. Thus we may express an arbitrary density operator as a linear combination of these 
operators as 

-1 d-l d-l 



P 



(1 



. b=0 m=l 1=0 



(5.6) 



For a given b, the operator XZ b possesses d eigenvectors, which we denote by \c: b) 



c = 0, 1, • ■ ■ , d — 1. In terms of the computational basis these eigenvectors are given by 



22| 



d-l 



b) = -= Va;^"- 1 )-^), XZ b \c-b) =u c \c]b). 
yd 



(5.7) 



n=0 



This equation defines d distinct bases (b = 0,1, ■ ■ ■ ,d — 1) which, when supplemented with 
the computational basis, which is an eigenbasis of the operator Z [see Eq. (15. lap ], forms a 
set of d + 1 MUB bases, i.e., 



(c;b\c';b) = 6 C , C ,, \{c; b\c'; b')\ 2 = 3, b±V, 
(n\n') = 8 ntn >, \(n\c;b)\ 2 



d' 
1 



d 

These equations can be proved straightforwardly by direct evaluation. 
We rewrite Eq. (15. 6p by adding and subtracting the m = terms as 

d-l d-l d-l 



(5.8a) 
(5.8b) 



P 



2 E E Tr [p((xz b rf] (xz b r - d 1 + J2 Tr [/W 

L b=0 m=0 1=0 



(5.9) 



The spectral representation of the operator XZ b is given by 



d-l 



XZ b = ^|c; b)uo c {c;b\. 



(5.10) 



c=0 



Note that the eigenvalues u c are non-degenerate. We obtain 



d-l 



Tr[p((XZ b ) m y] = Y,(c;b\p\c;b)u 



(5.11) 



c=0 
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TABLE I: Analogy between quantities for a discrete and a continuous Hilbcrt space 



Discrete case 


Continuous case 


X m Z l 




^ a , m(m — 1) , a a , 

X m Z z = u 2 b (XZ b ) m 


2 

giaigiftp g— i^-CS \^i(Cx+Sp) j r 


\c;b) 




XZ b \c;b) =uj c \c;b) 


gi(Gz+Sp)| x /.^ = g^'i^.^ 


(XZ b ) m \c;b) =uj mc \c-b) 





Hence the density operator of Eq. ( 15. 9ft takes the form 

d-l d-l 

p = 5^ I c ^)( c '^IpI c ^)( c ; & I + y^J n ) ( n \p\ n ) ( n \ -I - ( 5 - 12 ) 

6,c=0 n=0 

The matrix elements of p in the computational basis are then given by 

d-l 

(n'\p\n") = ^2 ( n '\^ & )( c ; 6 IpI c ; & )( c ; + (™'IpK)^>" - ( 5 -i3) 

6,c=0 

The density operator p is given in terms of probabilities, Eq. (15.121) . which are observable 
quantities; e.g., (c;b\p\c;b) is the probability to find the state \c;b) when the system is 
described by the density operator p. We thus find that p is reconstructed by using d + 



1 measurements [2l|. Each of these measurements yields d — 1 independent probability 
outcomes (since the probabilities add up to 1). This gives (d+ l)(d— 1) = d 2 — 1 quantities, 
which is precisely the number of independent parameters of a d- dimensional density matrix. 

Finally, we wish to call the reader's attention to the analogy between several quantities 
used in the present section and those introduced in the previous one, where a continuous, 
infinite-dimensional Hilbert space was used. This correspondence is indicated in Table [H 



VI. CONCLUSIONS AND REMARKS 



We have reviewed the approach to the quantum-state reconstruction problem based on 
the Wigner function and the Radon transform, pointing out its close analogy with classical 
tomography. We put emphasis on the role played by the Wigner function, which was shown 
to be analogous to that of the probability density in phase space for the classical problem. 
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The analysis underscores the intriguing fact that to reconstruct a quantum state we 
require the probabilities of all the phase-space plane, and not merely the probabilities along 
the position and momentum axes as might be implied by a positive reply to Pauli's query 
posed in the Introduction. 

Then we reviewed an alternative route for the state reconstruction which is based on 
the notion of mutually unbiased bases and does not make use of the Radon transform. We 
described its connection with the method based on the Wigner function. 

In addition, we showed that the concept of mutually unbiased bases can be applied to 
the problem of state reconstruction for a finite-dimensional Hilbert space, which is quite 
relevant for all applications to quantum computing. Finally, a parallel with the case of a 
continuous, infinite-dimensional Hilbert space is drawn. 
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Appendix A: Inverting the Radon Transform: proof of Eq. (12.31) 

Multiplying both sides of Eq. (12. 2p by e~ lkx ' and integrating over x' we find 



We identify the two sides of this equation with the Fourier transform pg(k) of po(x'), and 
the Fourier transform p(k x , k y ) of p(x, y), respectively, so that 




(Al) 





(A2) 



We recover p(x,y) as the inverse Fourier transform of p(k x , k y ), 




(A3) 
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where k x and k y are the Cartesian components of a wave number vector k; in polar coordi- 
nates we have 



k x = K cos <p, k y = K sin <fr, (A4a) 
K — |k| > 0. (A4b) 

The density p(x, y) becomes 

I roc p2tt 

p(x,y) = T - — - dK K d<Pe iK{xcos,p+ysin ^p(Kcos4),Ksm<i)). (A5) 
(2vr) 2 Jo Jo 

While the variable k in Eq. (1A2|) is defined in the interval (—00, 00), the radial variable K in 
Eq. (1A5I) is defined to be non-negative and in the interval (0, 00). The range of integration 
of the variable K can be extended to the full real axis by first splitting the interval of 
integration of 4> into the intervals (0, n) and (tt, 2tt) and then making the change of variables 
= (j)' + 7r, K = —K' in the integral over the second interval, to obtain 

1 /»00 PIT 

p( x ,y) = ^^/ dk\k\ d6 e iHxcos9+ysin0) p(k cos 6, k sin 6). (A6) 
{2%y J^ J 

We identify the last factor with the quantity p~e(k), Eq. ( 1A2I) . and substitute pe(k) from the 
left-hand side of Eq. ( 1A1I) to write 



I poo pn /*oo 

p( x ,y) = ^^/ dk\k\ d Q e ik(* cos e+y sine) / dx' e~ ikX ' P g(x') . (A7) 

(271T J.oo Jo J-oo 

Defining the integral 

/oo 
|fc|e- ifc «-l fc l^, 6>0, (A8) 
-00 

and identifying £ = a?' — (xcos^ + ysin#), we write Eq. (IA7I) as 

p(x,y) = lim d8 dx' fJx' - (xcos6 + ysinfl)) pgix'). (A9) 

Thus our task is to study the function / € (£), which we write as 

/0 ;>oo 
(-A;)e- ifcK+ie) ^+ / ke~ ik ^- ie) dk 
00 Jo 

0/ 1 1 \ flam (A10) 



where 



^(0 = 2^^. (All) 
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Using the abbreviation a = x cos 6 + y sin 6, we write the last integral in Eq. (1A9|) as 

pe(x')f e (x' - a)da/ = - f°° dPe } X \ e {x' - a)dx' (A12) 



dx' 



where we have used the definition (1A10I) and we have integrated by parts, assuming the 
integrated term to vanish for sufficiently large values of the argument. 

The function g e {£) is shown schematically in Fig. [3] As e — > 0, the integral of Eq. ( 1A12I) 







-£ 









FIG. 3: Schematic plot of the function g e (£) defined in Eq. (jAlljl 
tends to the principal-value integral 



limJ 6 = -2V 



dp e (x')/dx' ^ 



x' — a 



Substituting this result in Eq. flA9j) . we then find Eq. (I2.3P in the text. 



(A13) 



Appendix B: Proof of Eq. (I37L4D 

We first remark that it is easy to prove the operator identity 

P£, = 5{x' - Xg) . 



(Bl) 



Therefore, we compute the required Wigner transform of the projection operator (IBlj) as 



Q + - 
q 2 



S[x' - {Cq + Sp)} 



q -2 )e 



-tpy 



dy 



(B2a) 
(B2b) 
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It is convenient to work with the Fourier transform of this last expression with respect to 
the variable x'\ i.e., 



e ikx W r e / (q,p)dx' 

ik(Cq+Sp) 



q+ 2 



-ipfj 



dy 



g ~2 k^ 1 SC 



q + 



y 



^ikCq^ikSp 



y 



-ipy 



dy 



(B3a) 
(B3b) 
(B3c) 



where use was made of the Baker- Campbell-Hausdorff (BCH) identity (Ref. p. 442) 

e A+B = e A e B e -l[A,B] ^ (B4) 

valid for any two Hermitean operators A, B, whose commutator commutes with each of 
them, i.e., [A, [A,B]] = [B, [A,B]] = 0. Introducing inside the matrix element of Eq. (1B3cP 
a complete set of eigenstates of position and of momentum right after the first and second 
exponentials, respectively, we find 

e ikx 'W ¥ e (q,p)dx' = e ik{Cq+Sp) . (B5) 

x 1 

The inverse Fourier transform of this last expression gives the result of Eq. f )3.14p . 



Appendix C: Proof of Eq. (14.71) 

To prove Eq. (14.71) we impose the requirements of Eq. (14. 6 p on the solution, Eq. (14. 5p . 
1) The ortho- normalization condition, Eq. (I4.6a[) . imposed on the wave function ip x ' t e(x) 
of Eq. (14.51) gives, for the function F(x', 9), 

p i<t>e{x') 

F&,0) = -=—■, (CI) 
y 2tt \ smv\ 

4>e{%') is an arbitrary phase, dependent on x' and 9. The wave function ip x ' t e(x) becomes 

w, , J x ) = 1 p ~^[(x^co S e-2xx') +l M^)} ( C2 ) 
^tt I sine I 



2) The first requirement in Eq. ( 14.6bl) for the matrix elements of Xg with the wave function 
of Eq. (1C2I) is automatically fulfilled, since i) our starting point has been the eigenvalue 
equation, Eq. (13. lip , and ii) the wave function of Eq. ( 1C2|) satisfies the orthonormalization 
condition, Eq. (14. 6 aft . 
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3) We first compute the matrix element of Pg which appears on the left-hand side of the 
second requirement in Eq. (14.6bj) . Using the definition of Pg given in Eq. (14. 2p and Eq. 
( IC2j) we have, in the coordinate representation 



(xi,6\P e \x 2 , 



d 

sin 9 x — % cos 9— 

ox 



i) X2j e(x)dx, 



cos 9 



:X 



sin 9 sin 9 



x 2 



i) X2) e{x)dx, 
cos 9 



-i e 4 ^(^)-0«(-i)] ( 5'( Xl _ X2 \ + * X2 S ( Xl _ X2 \ 

sin 9 



(C3) 



We write, for the above exponential, the Taylor expansion 

e i[(j) B {x2)-<j> B {x 1 )} _ y _|_ i 



x 2 - xM'Jxt) H 4>q(xx 



1 

~2! 
+ ••• , 



(x 2 - xxYW^xt))* + 2 



2! 

(x 2 - Xi) 3 



2! 



^(xi)^(xi) + 



(C4) 



where the primes mean derivatives with respect to the argument. We use the 5-function 
identities xS'(x) = —S(x), x n 5'(x) = 0, n > 2, to write the matrix element, Eq. ( 1C3I) . as 



cos 9 



-i6'(x! - x 2 ) + cf) g {xi)8(xi - x 2 )] + -t— qX\5(xx - x 2 ). 

sin u 



(xi,9\P e \x 2 ,9) 

In order to satisfy the second requirement in Eq. (14.6bj) we thus need 
with the solution 



cos 9 , 

x 



sin 9 

cos 9 x' 2 



+ ¥>(*), 



sin 9 2 

where <p(9) is an arbitrary function of 9. 

The wave function ip x \e{x) of Eq. (1C2I) then becomes 

e i<p(6) 



(C5) 



(C6) 



(C7) 



4) Choosing, for the phase 



27T sin^ I 



7T 



e -OT[( a:2 + a: ' 2 ) cos0 - 2a:a: '] 



(a 



<P(0) = ^sgn(sini 



9 

2' 



(C9) 



we satisfy the requirements of Eq. (I4.6cl) . 

We finally find the wave function of Eq. (14. 7p 
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Appendix D: The matrix elements of the operator exp(i9h) 



In this Appendix we compute the matrix elements of the operator exp (i8n) in the original 
basis \x). We have 



x 2 +x ,2 e iy ) n 



(Dla) 
(Dlb) 



where ip n {x) are the one- dimensional harmonic oscillator wave functions (see Ref. 27J, p. 
61, Eq. (5.24), where the variable x has been replaced by the dimensionless xj J as used 
in this paper) 



-e x 2 H n (x), 



(D2) 



7r 1 / 4 V2 n n\ 

H n (x) being Hermite polynomials. We compute the sum in Eq. ( IDlbj) using the identity 
(Ref. [28|, p. 781, Problem 6.12) 



t n 1 

EH n (x)H n (x') = exp 

r> v 



2xx't — t 2 (x 2 + x 
1-t 2 



(D3) 



with the result 



{x\e l0n \x') = e 



e 2 



[(x 2 +x' 2 ) cos 6-2xx'} 



^2vr(-z) sin^ 

This result is identical to that of Eq. ( 14.91) if we choose, for the square root, the branch 



(D4) 



sin 6 



e l 4 a/] sin 9\, for sin# > 
e +I 5^/| sin#|, for sin 6* < 



"^ sgn(sine) V^hilj. 



(D5) 
(D6) 



Appendix E: A simple way to derive the result (14.14|) 

It will suffice to evaluate the quantity \(x\x';9)\ 2 ; this is the probability to find x in a 
unit interval around the value x when the system has been prepared in the state \x'; 9). We 
find 

|(x|x';#>| 2 = (x\F x '-,e\x) = Tr{F x -eh) (Ela) 
= // W K ^q,p)W tm (q,p)^ , (Elb) 
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where in the last line we have used Eq. ( 13. 7p to express our probability in terms of Wigner 
transforms. The Wigner transform of the projector F x '-e is found from Eq. (I3.14p and that 
for the projector F x is simply 5(x — q). We thus write the last equation as 

\{x\x'-6)\ 2 = J j 5(x' -(qC + pS))5(x-q)^ (E2a) 

= — / S(x' - (xC + P S))dp . (E2b) 
2n J 

On the one hand, this integral can be evaluated directly, giving 

ww' = £f&{*-^)*=m- (E3) 

just as in Eq. ( 14. 14ft with 9\ = 9 and 62 = 0. On the other hand, the appearance of the 
factor \S\ in the denominator of the result (1E3I) can be understood by using an intuitive 
geometrical argument starting from (1E2bl) . as follows. We approximate the delta function 
occurring in Eq. (lE2bl) by the step 

f ± if x' e (x' - 5 -f,x' + 5 -f) 
5(x — (xC + pS)) « usx'(x ) = < 1 , (E4) 

[ 0, if x' i (x> - s -f,x' + 5 -f) 

the delta function being attained in the limit 8x' — > 0. The non-zero region is indicated as 
the shaded area in Fig. HJ The segment along the p-axis, over which we are integrating, 
contained inside that area, is Sp = 5x' / sin 9. The integral in Eq. (lE2bj) is thus given by 

1 8x' 1 



5x' sin 9 sin 9 
thus reproducing the formal result of Eq. ( 1E3D . 



(E5) 
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